Your instructors
Michiel Stock
UGent postdoc computational biology
Authors and contributor several Julia
packages, founder of Julia Doctoral
Schools
Maxime Van Haeverbeke
UGent PhD student data science for
bioimpedance
Author of Julia packages equivalent
electrical circuits
Goals
1. Learning what Julia is and why you should use it
2. Setting up Julia and learn to work in the Pluto notebook environment
3. Define, solve and fit ordinary differential equations in Julia
4. Perform Bayesian inference using probabilistic programming
5. Give you pointers to start your own Julia journey!
Oddish will indicate an exercise
We want a language that’s open source, with a liberal license. We want
the speed of C with the dynamism of Ruby. We want a language that’s
homoiconic, with true macros like Lisp, but with obvious, familiar
mathematical notation like Matlab. We want something as usable for
general programming as Python, as easy for statistics as R, as
natural for string processing as Perl, as powerful for linear algebra as
Matlab, as good at gluing programs together as the shell. Something
that is dirt simple to learn, yet keeps the most serious hackers happy.
We want it interactive and we want it compiled.
(Did we mention it should be as fast as C?)
Jeff Bezanson, Stefan Karpinski,
Viral B. Shah, Alan Edelman
The secret sauce of
JuliaLang
Type-based dispatch allows for the generation of
highly optimised and specialised methods
(functions) for a given input combination.
Julia gives you super powers!
Roesch, E., Greener, J. G., MacLean, A. L., Nassar, H., Rackauckas, C., Holy, T. E., & Stumpf, M. P. H. (2021). Julia for Biologists. Retrieved from http://arxiv.org/abs/2109.09973
1. Speed
2. Abstraction
3. Metaprogramming
Three main strengths of Julia:
1. Julia is fast!
●Typically much faster than Python, R, …
●Before a function is run, it is compiled
from the given type signature
●Close to the speed of C, if code is
properly optimized (not necessarily easy)
●Can be made faster than C and the like,
because one can use high-level tricks
such as metaprogramming
Julia is fast!
Comparison of gradient descent in Python and Julia
Algorithm for solving
Comparison: gradient descent in Python and Julia
Python Julia
Comparison: gradient descent in Python and Julia
Julia
Python
2. Abstraction by type system and dispatch
Julia allows users to construct their own types
within the native type system.
Multiple dispatch allows to specify certain
functions for certain combinations of types.
This is extremely powerful!
Composability
Code of different sources can easily be combined into novel ways:
- Zygote.jl + Colors.jl = differentiation on colors
- DifferentialEquations.jl + Measurements.jl = ODEs with error bars
- DifferentialEquations.jl + Flux.jl = neural ODEs
Example of abstractions: Measurements.jl
julia> using Measurements
julia> a = 4.5 ± 0.2
4.5 ± 0.2
julia> b = 3.8 ± 0.4
3.8 ± 0.4
julia> a + b
8.3 ± 0.45
julia> a * b
17.1 ±2.0
julia> cos(a)
-0.21 ±0.2
julia> a / a
1.0 ± 0.0
3. Metaprogramming
Metaprogramming is using a programming language to write new code.
Julia is homoiconic: Julia code can be represented as a Julia data structure
(similar to languages such as LISP)
Pluto notebooks
DifferentialEquations.jl
Installing Pluto notebooks
1. Press the “]” key in julia to get into the package manager interface
2. Install the package by typing “add Pluto”
3. Import the module by typing “using Pluto”
4. Open the notebooks interface by typing the command “Pluto.run()”
Exercise: getting started with the Pluto notebooks
Start by getting familiar with Julia using the “Basic mathematics”
sample notebook.
At the same time, load the other notebooks of the workshop. These
can be loaded by first Ctrl+clicking on the Pluto.jl icon at the top
right of the screen and navigating to the 1_ODEs.jl and 2_catalyst.jl
notebooks.
See you back at 16H00 after the break!
Ordinary Differential Equations
DifferentialEquations.jl
Modelling with ordinary differential equations
Most (interesting) ODEs can only be solved numerically!
Differential equations can be derived from conservation laws
accumulation
in out
system
Change through time = input - output + production - degradation
Example: logistic growth
Logistic growth models the growth of a population with a
fixed capacity.
Growth
rate
Carrying
capacity
Solve the ODE using DifferentialEquations.jl
Note: works with all kinds of solvers for
when your problem is stiff etc.
Plot the solution
Exercise: a variable growing rate
Consider a slightly more different model
Where the growth fluctuates during the day. There is an extra parameter wthat
controls the periodicity of the growth. Implement this ODE as
g(x, (μ, K, w), t)
Plot the result and play with w.
Calibrate the parameters to data
How can we choose our
parameters so that our model
nicely follows these
observations?
Couples of (time, yield)
Calibrating the model
Consider the sum of squared residuals between model
and data
Find the parameters by minimizing this error
ODEs with multiple variables
Consider the R* model with multiple species that share a
common resource
Crop and weeds, competing for nitrogen
Solution of the R* model
Modify the R* model
Reaction systems
Catalyst.jl
Domain-specific languages (Catalyst.jl)
How to simulate a simple reaction network
Modelling photosynthesis
Incorporating photorespiration
Going from C3 to C4 models: spatial compartiments
Implement your own reaction
Probabilistic modeling
Turing.jl
Bayesian reasoning
Schoot, R., Depaoli, S., King, R., Kramer, B.,
Mamp, K., Tadesse, M. G., … Yau, C.
(2021). Bayesian statistics and
modelling. Nature Reviews Methods
Primers, 1(1), 1.
posterior
evidence
likelihood prior
Computing the evidence is usually the hard part!
Approximating the posterior with samples
Exploring the posterior space by clever random
walks, generating samples from the posterior.
Almost any inference question can be answered by computing
expected values based on posterior samples!
Probabilistic programming
Using a programming language to simulate samples from
arbitrary complex probability distributions.
Probabilistic programming languages
Example: how effective is a treatment?
We are testing two potential treatments for
flower induction. Treatment A is effective for
29 in 30 plants, treatment B is effective for 17
in 20 plants. Which treatment is more
effective?
@model function treatment(N_A₊, N_A₋, N_B₊, N_B₋)
# prior distributions of probability of success
p_A ~ Uniform(0, 1) # probability that A is effective
p_B ~ Uniform(0, 1) # probability that B is effective
# number of cases for each treatment
N_A = N_A₊+ N_A₋
N_B = N_B₊+ N_B₋
# sample successes
N_A₊~ Binomial(N_A, p_A)
N_B₊~ Binomial(N_B, p_B)
end
julia> p_A_better_than_B = mean(chain[:p_A] .> chain[:p_B])
0.931
julia> p_A_much_better_than_B = mean(chain[:p_A] .> 1.2chain[:p_B])
0.327
Bayesian curve fitting
Recall our data of dry mass through time:
Given two growth models (3 parameters
each):
- Gompertz
- Logistic
Comparing the chains of the two growth models
Gompertz Logistic
Comparing the chains of the two growth models
Gompertz Logistic
Comparing the chains of the two growth models
Gompertz Logistic
Sampled parameters for the logistic model
Further resources
